† Corresponding author. E-mail:
Transcranial focused ultrasound is a booming noninvasive therapy for brain stimuli. The Kelvin–Voigt equations are employed to calculate the sound field created by focusing a 256-element planar phased array through a monkey skull with the time-reversal method. Mode conversions between compressional and shear waves exist in the skull. Therefore, the wave field separation method is introduced to calculate the contributions of the two waves to the acoustic intensity and the heat source, respectively. The Pennes equation is used to depict the temperature field induced by ultrasound. Five computational models with the same incident angle of 0° and different distances from the focus for the skull and three computational models at different incident angles and the same distance from the focus for the skull are studied. Numerical results indicate that for all computational models, the acoustic intensity at the focus with mode conversions is 12.05% less than that without mode conversions on average. For the temperature rise, this percentage is 12.02%. Besides, an underestimation of both the acoustic intensity and the temperature rise in the skull tends to occur if mode conversions are ignored. However, if the incident angle exceeds 30°, the rules of the over- and under-estimation may be reversed. Moreover, shear waves contribute 20.54% of the acoustic intensity and 20.74% of the temperature rise in the skull on average for all computational models. The percentage of the temperature rise in the skull from shear waves declines with the increase of the duration of the ultrasound.
Transcranial focused ultrasound (tcFUS) is a rapidly developing noninvasive technique for brain diseases, with application prospects in noninvasive thermal ablations of intracranial tumors,[1] neuromodulations,[2,3] treatment of neurological disorders,[4,5] blood brain barrier opening,[6,7] intracranial targeted drug delivery,[4,8] and intracranial thrombolysis.[9,10] The skull on the wave propagation path will easily cause the overheating of itself and a focal shift of ultrasound, as the skull has much larger impedance difference, much heavier absorption, and much stronger inhomogeneity than its surrounding media.[11–13] Thus, it is necessary to understand and then numerically reproduce the propagation of tcFUS and the temperature rise effect induced by tcFUS, for the purpose of better guiding the tcFUS applications.
In the early stage before 2000, the tcFUS studies only considered compressional waves while neglecting shear waves, on the assumption that shear waves did not play a significant role in tcFUS, especially when the incident angle of ultrasound was smaller than 20°.[14–16] The linear and nonlinear acoustic wave equations in fluids, commonly used in high-intensity focused ultrasound (HIFU), were adopted to reproduce the tcFUS field.[17–22] However, since 2000, mode conversions between compressional and shear waves in the skull have attracted increasing attention and research in tcFUS. By utilizing a simple layered model containing the skull, Clement et al.,[15,23] and Hayner and Hynynen,[16] pointed out that mode conversions occur in the skull and shear waves play an important role, especially when the incident angle is greater than 20°. White et al. measured that in the skull, the shear wave velocity is about half of the compressional wave velocity and the attenuation coefficient of shear waves is much higher than that of compressional waves on average.[24] Pinton et al. further showed that the absorption coefficient of shear waves is greater than that of compressional waves in the skull and that the acoustic attenuation caused by the acoustic absorption is less than those caused by the combination of reflections, scatterings and mode conversions.[25] Pichardo and Hynynen employed the Rayleigh–Somerfeld integral for a multilayer case to numerically study the treatment of near-skull brain tissue using shear-mode conversion.[26] Pulkkinen et al. investigated the skull base heating and the tcFUS therapy for the chronic neuropathic pain based on the coupled fluid-solid wave equations and the Pennes bio-heat equation, and the numerical results showed good agreement with the experimental results.[27,28] Song et al. studied standing-wave formation in a human skull with the coupled fluid–solid wave equations.[29] In recent years, the wave equations in solids such as the Kelvin–Voigt equation and the Biot equation have been employed to reproduce the tcFUS field, which have achieved better agreement with the experimental results than the theoretical results based on the aforementioned wave equaitons in fluids.[30–35]
Previous work including mode conversions in the skull mainly focused on the overall distribution of the tcFUS field generated by a spherical phased array.[26–29,31,34] Pulkkinen et al. further studied the temperature field caused by the overall tcFUS field.[27,28] These researches lacked the exploration of the separated contributions from compressional waves or shear waves to the tcFUS and temperature fields, especially in the skull. However, in this paper, we will numerically focus a 256-element planar phased array through a monkey skull and introduce the wave field separation method to the Kelvin–Voigt equation and the Pennes equation. The influences of mode conversions in the skull on the tcFUS and temperature fields are discussed, the contributions from compressional and shear waves to the acoustic intensity and the temperature rise in the skull are also analysed, respectively.
The Kelvin–Voigt equation in isotropic viscoelastic solids has been successfully employed to simulate the tcFUS field with mode conversions in the skull.[32,34] In an elastic material, the strain at each point is only dependent on the instantaneous local stress. The stress and strain are related by the stiffness. For an anisotropic medium, this relationship can be written using Einstein summation notation as
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn1.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn2.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn3.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn4.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn5.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn6.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn7.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn8.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn9.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn10.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn11.gif)
![]() |
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn12.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn13.gif)
To describe the propagation of elastic waves, the momentum conservation relationship between stress and strain must be combined, which is written by the stress and the particle velocity
![]() |
Equation (
![]() |
Equations (
The wave field separation is based on the fact that a compressional wave has no vortex and a shear wave has no divergence either. The particle velocity is separated into contributions of compressional and shear waves denoted respectively as
![]() |
![]() |
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn17.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn18.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn19.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn20.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn21.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn22.gif)
![]() |
Since the compressional and shear parameters of the skull are different, the heat deposition in the skull is contributed by the separated compressional and shear waves above. For single-frequency ultrasound used in tcFUS, the heat deposition can be written as
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn23.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn24.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn25.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn26.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn27.gif)
The Pennes bio-heat equation[37] is commonly used to describe the temperature field and can be written as
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn28.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn29.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn30.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn31.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn32.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn33.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn34.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn35.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn36.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn37.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn38.gif)
The staggered-grid finite difference time domain (FDTD) method[38,39] is employed to numerically solve the Kelvin–Voigt equation and the wave field separation equations. The parameter averaging technique[40] is used to satisfy the continuity conditions on the interfaces of solids (skull) and fluids (water and brain inside the skull). The non-splitting convolutional perfectly matched layer (NCPML)[41,42] is utilized to eliminate the numerical reflections at the numerical boundary. The FDTD method is also used to numerically solve the Pennes equation.
Generally, computed tomography (CT) scans of a monkey skull are used to reconstruct the three-dimensional (3D) skull model, as shown in Fig.
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn39.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn40.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn41.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn42.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn43.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn44.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn45.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn46.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn47.gif)
The planar phased array for tcFUS is shown in Fig.
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn49.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn50.gif)
![]() | Fig. 2. (color online) Schematic diagrams of (a) planar phased array, and (b) cross section across central axis with focus. |
Eight computational models are chosen as shown in Fig.
![]() | Fig. 3. (color online) Computational models where ((a)–(e)) θ = 0° and FD increases from 10 mm to 50 mm in steps of 10 mm, and ((f)–(h)) FD = 30 mm and θ increases from 10° to 30° in steps of 10°. |
The time-reversal method[44] is used to precisely focus ultrasound from the array through the skull at the focus. Ultrasound emitted from a virtual sinusoidal point source
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn53.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn54.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn55.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn56.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn57.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn58.gif)
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn59.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn60.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn61.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn62.gif)
The initial temperature is
![]() | Table 1.
Acoustic and thermal parameters used in the simulation. . |
What needs to be pointed out is that complex guided waves may be excited in the skull, because the thickness of the skull is close to the wave length, multiple reflections and refractions occur at the boundaries of the skull and the inhomogeneity of the skull leads to strong mode conversions and coupling between compressional and shear waves in the skull. Nevertheless, the guided wave is another form of compressional wave or shear wave in the skull, which has been included in the Kelvin–Voigt equation as the Kelvin–Voigt equation can describe the full wave field in isotropic solids. The guided waves in the skull may make sense in monitoring the tcFUS therapy and will not be discussed here since what we are going to study in this paper is not related with the guided waves.
The normalized acoustic intensity (
From Figs.
According to Figs.
Figures
Figure
For the 256-element phased array with a uniform emitting acoustic pressure amplitude, numerical results of the normalized acoustic intensities above indicate that a focus designed to be closer to the skull or a larger central incident angle from the array to the skull make it easier to cause serious acoustic energy to deposit around the skull, which is detrimental to the tcFUS therapy. For a central incident angle smaller than 30°, an overestimation of acoustic energy deposition around the focus and an underestimation of acoustic energy deposition around the surfaces of the skull will simultaneously occur when mode conversions are neglected. Nevertheless, when the central incident angle exceeds 30°, the situation may be reversed. We suggest that a location deep in the brain combined with a central incident angle close to 0° may be more suitable for ultrasound emitted from the array to be focused at.
The duration of the tcFUS field is set to be TD = 1 s and the emitting pressure is chosen to be
Like SubSection 3.1, the temperature rise at the focus and the maximal temperature rise around the skull for each of all computational models are picked up from Figs.
First of all, figures
Figure
In Fig.
Likewise, an overestimation of the temperature rise around the focus and an underestimation of the temperature rise around the surfaces of the skull simultaneously exist for a central incident angle smaller than 30° when mode conversions are neglected. A reversed situation may also occur if the central incident angle is beyond 30°. Nevertheless, comparing with numerical results of the acoustic intensities in SubSection 3.1, the big characteristic or difference for numerical results of the ultrasound-induced temperature rises is that the heat deposition problem is greatly exacerbated due to the strong absorption of the skull to the acoustic energy. Therefore, such a planar phased array with a uniform emitting acoustic pressure amplitude is more applicable to the non-thermal tcFUS therapy in order to avoid the skull heating problem, especially when the focus is deep in the brain and the central axis of the array is perpendicular to the skull as much as possible. As for the thermal tcFUS therapy utilizing this array, it is better that the microbubbles should be used to deliberately enhance the energy deposition at the focus, in order to relieve the skull heating problem.
Equation (
From Figs.
Moreover, figure
![]() |
![](cpb_27_2_024302/cpb_27_2_024302_ieqn94.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn95.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn96.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn97.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn98.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn99.gif)
![](cpb_27_2_024302/cpb_27_2_024302_ieqn100.gif)
According to Fig.
Furthermore, as just mentioned above, the difference between the average proportions of the acoustic intensity and the temperature rise in the skull from shear waves is indeed existent yet tiny and inconspicuous, as shown in Fig.
Taking the computational model of FD = 30 mm and θ = 0° for example, we give the spatial distributions of the proportions of the temperature rises induced by shear waves with TD = 0.01, 0.5, 2.0, and 3.0 s in Fig.
In this work, to account for the different absorptions of the skull to compressional and shear waves, the wave field separation method is introduced into the Kelvin–Voigt equation and the Pennes equation for the study of the influences of mode conversions in the skull on the tcFUS and temperature fields.
From the simulation results of focusing a 256-element planar phased array through a monkey skull with eight computational models, it is confirmed that mode conversions in the skull do play an indispensable role in the tcFUS and temperature fields. Shear waves induce considerable acoustic intensity and temperature rise in the skull. Ignoring mode conversions, i.e., ignoring shear waves, will lead to a poor assessment of the tcFUS field and the tcFUS-induced temperature field, which will be harmful to guiding the tcFUS applications. In addition, a focus close to the skull or a large incident angle from the array to the skull is apt to induce the overheating of the skull, which is against the noninvasive aim of the tcFUS therapy. Focusing the ultrasound at a location deep in the brain with a small incident angle is better.
By considering the mode conversions in the skull and separately calculating the contributions from compressional and shear waves to the skull with different focusing models, this work can provide a better evaluation and a better guidance for the practical experiments or treatments by utilizing the tcFUS.
The parameters of the 256-element planar phased array and the CT scans of the monkey skull were provided by Shenzhen Institutes of Advanced Technology, Chinese Academy of Sciences. The numerical calculations were performed on the “YUAN” supercomputing clusters of the Supercomputing Center of the Chinese Academy of Sciences.
[1] | |
[2] | |
[3] | |
[4] | |
[5] | |
[6] | |
[7] | |
[8] | |
[9] | |
[10] | |
[11] | |
[12] | |
[13] | |
[14] | |
[15] | |
[16] | |
[17] | |
[18] | |
[19] | |
[20] | |
[21] | |
[22] | |
[23] | |
[24] | |
[25] | |
[26] | |
[27] | |
[28] | |
[29] | |
[30] | |
[31] | |
[32] | |
[33] | |
[34] | |
[35] | |
[36] | |
[37] | |
[38] | |
[39] | |
[40] | |
[41] | |
[42] | |
[43] | |
[44] | |
[45] |